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We investigate the ground and low excited states of a rotating, weakly interacting Bose-Einstein 
condensed gas in a harmonic trap for a given angular momentum. Analytical results in various 
limits, as well as numerical results are presented, and these are compared with those of previous 
studies. Within the mean-field approximation and for repulsive interaction between the atoms, we 
find that for very low values of the total angular momentum per particle, L/N — » 0, where Lfi is 
the angular momentum and TV is the total number of particles, the angular momentum is carried by 
quadrupolar (\m\ = 2) surface modes. For L/N = 1 a vortex-like state is formed and all the atoms 
occupy the m — 1 state. For small negative values of L/N — 1 the states with m = and m = 2 
become populated, and for small positive values of L/N — 1 atoms in the states with m = 5 and 
m — 6 carry the additional angular momentum. In the whole region < L/N < 1 we have strong 
analytic and numerical evidence that the interaction energy drops linearly as a function of L/N. We 
have also found that an array of singly quantized vortices is formed as L/N increases. Finally we 
have gone beyond the mean-field approximation and have calculated the energy of the lowest state 
up to order N for small negative values of L/N — 1, as well as the energy of the low- lying excited 
states. 

PACS numbers: 03.75.Fi, 05.30.Jp, 67.40.Db, 67.40.Vs 



I. INTRODUCTION 

One of the basic questions about Bose-Einstein con- 
densates in trapped alkali atom vapors |j| is how they 
behave under rotation. A lot of theoretical work has 
been done on this subject, both analytical, and numer- 
ical and the problem has been studied theoreti- 
cally, both in the Thomas-Fermi limit of strong interac- 
tions |||]|§ and in the limit of weak interactions , 
which we consider in this paper. 

In Ref. H Wilkin et al. considered a weakly interact- 
ing Bose gas with attractive interactions and showed that 
in the lowest energy state of a given angular momentum, 
the angular momentum is carried by the center of mass 
motion. Butts and Rokhsar calculated numerically the 
moment of inertia of a weakly-interactingtrapped Bose 
gas with effective repulsive interactions m. One of us 
identified the elementary modes of excitation for small 
angular momentum and demonstrated in Ref. || that 
a system of rotating weakly-interacting bosons exhibits 
two additional kinds of condensation associated with the 
nature of low-lying excitations. Finally Bertsch and Pa- 
penbrock performed in Ref. Q] exact numerical diagonal- 
ization within the subspace of states with a given angular 
momentum, which are degenerate in the absence of in- 
teractions. 

Experimentally the detection of vortex states in a two- 
component system has been reported by Matthews et al. 
||, while Madison et al. have provided evidence for 
the formation of vortex states in a stirred one-component 
Bose-Einstein condensate. 

Our basic goal in this study is to identify the lowest 
energy states of a harmonically trapped, weakly inter- 
acting Bose gas for a given angular momentum L. As 



we show below these states are selected by the interac- 
tions. In Sec. II we describe the model and discuss the 
degeneracy of the many-body states for a given angular 
momentum in the absence of interactions. In Sec. Ill we 
use the mean-field approximation to calculate the inter- 
action energy, and derive numerical and analytical results 
under various conditions. In Sec. IV we describe how one 
can go beyond the mean-field approximation and study 
as an example the specific case of small negative L/N—l. 
Finally in Sec. V we give our conclusions. 

II. THE MODEL 

Our starting point is the Hamiltonian H, given by 

H = H + V. (1) 

Here 

H ° = E -|^ v * + E \ M ^ + yf) + /(*) ( 2 ) 

i i 

includes the kinetic energy of the particles and their po- 
tential energy due to the trapping potential. The axis 
of rotation is taken to be the z axis, and the trapping 
potential is assumed to be that of an isotropic harmonic 
oscillator of frequency lo in the x-y plane. Also M is the 
mass of the atoms. Our results do not depend on the 
trapping potential f(z) in the z direction. The interac- 
tion V between atoms is assumed to be of zero range, 

V=lu J2S(ri-rj), (3) 

#3 
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where Uq = Airh 2 a/M is the strength of the effective two- 
body interaction, with a being the scattering length for 
atom-atom collisions. We assume that the interaction 
is repulsive, a > 0. Attractive interactions have been 
studied in Refs. H|). 

Much theoretical work on rotating condensates has 
been done in the Thomas-Fermi limit of strong interac- 
tions, where the superfluid coherence length 



-1/2 



(4) 



n being the particle density, is much less than the size of 
the cloud. Under these conditions the system is expected 
to exhibit superfluid properties much like those of liquid 
helium II |11|. In this study we examine the opposite 
limit of weak interactions, nllo <C twj and nUo <C AE Z , 
where A.E Z is the energy separation between the first 
excited state and the ground state for motion in the z 
direction. Under the above conditions 



aw V Na J 



(5) 



where N is the number of atoms in the trap, a osc — 
(Ti/Moj) 1 / 2 is the oscillator length, and a z is the charac- 
teristic length associated with the motion of the atoms 
along the z axis. Therefore the coherence length is larger 
than the size of the cloud, and the situation is analogous 
to that for BCS pairing of nucleons in nuclei. 

Since we consider rotation around the z axis, the con- 
dition nUo <C AE Z implies that the motion along this 
axis is frozen out and the problem is essentially two- 
dimensional. It is well known that for the harmonic 
oscillator potential in two dimensions the single-particle 
energies e are given in the absence of interactions by 



(2n r 



(6) 



where n r is the radial quantum number, and m is the 
quantum number corresponding to the angular momen- 
tum. In the lowest energy state of the many boson system 
all particles are in states with n r = 0, and with m being 
zero or having the same sign as the total angular mo- 
mentum. The energy of the lowest state of a system of 
non-interacting bosons with angular momentum L mea- 
sured relative to that of the ground state is therefore 



E = Lhu. 



(7) 



There is a huge degeneracy corresponding to the many 
different ways of distributing L quanta of angular mo- 
mentum among N atoms. Interactions between the 
atoms lift the degeneracy. We incorporate the effect of 
the interactions in both the mean-field approximation, 
as well as by diagonalization within some appropriately 
chosen truncated space of degenerate states. We describe 
the two methods separately below. 
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FIG. 1. The numerical result for |c m | which comes from 
the minimization of the energy, as a function of L/N. The 
numbers refer to the corresponding states with angular mo- 
mentum mh. The lowest two curves in the region L/N ~ 1.5 
give the occupancy of the m = 5 (higher) and m = 6 (lower) 
states. 




FIG. 2. The interaction energy, (V), in units of TV 2 wo as a 
function of L/N. 



III. MEAN-FIELD APPROXIMATION 

We start with the mean-field Gross-Pitaevskii ap- 
proach. Butts and Rokhsar have used this method to 
derive numerical results for the moment of inertia of a 
Bose gas Q]. In this scheme the many-body condensate 
wavefunction with N particles and L units of angular mo- 
mentum ^ L../v( r ij r 2j ■ • ■ i ?n) is taken to be the product 
of the single-particle states ^(r,), 



L.N 



(ri, r 2 , . . . , rjv) = *(ri) x *(r 2 ) . . . *(rjy). (8) 



The single-particle states ^(r,) can be expanded in terms 
of the harmonic-oscillator eigenstates § m (ri): 



2 



(9) 



m=0 



where the c m are variational parameters, which are com- 
plex in general and are functions of L. The summation in 
Eq. (|^) is restricted to positive to, since states with neg- 
ative m do not belong to the space of degenerate states. 
The quantity |c m | 2 gives the occupation probability for 
state $„, . Also 



$ m (r) 
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(TOlTraL) 1 / 2 
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(10) 



Here p, z, and (/> are cylindrical polar coordinates. In the 
above expression we have assumed that the bosons are in 
their ground state g(z) along the axis of rotation. The 
expectation value of the interaction energy V in the state 
given by Eq. (@) is 



(V) = ~N(N-1)U 



l*| 4 dr. 



(11) 



To find the lowest energy state we calculate (V) as a func- 
tion of the variational parameters c m , and minimize it 
with respect to them under the following two constraints: 
the normalization condition, 



i. 



(12) 



and the condition that the expectation value of the an- 
gular momentum per particle be fixed, 



L/N. 



(13) 



The parameters c m are complex in general, and there- 
fore both their magnitudes, and their phases need to be 
determined. However Eqs. ([l2]) and (|l3| ) impose two con- 
straints on the magnitudes of the c m . Furthermore, the 
overall phase of the wavefunction is arbitrary. Finally the 
rotational symmetry of the confining potential implies 
that the origin of the angular coordinate is arbitrary, 
which corresponds to the condition for conservation of 
angular momentum, which holds even in the presence of 
interactions. Therefore if the expansion ^ is truncated 
at a value TO max , the number of independent variables is 
2 x (m max + 1) - 4 = 2(m max - 1). 



A. Numerical results 

We have examined the problem numerically with TO max 
up to 9. The total number of terms in the expression for 
(V) is 125 in this case. The result of such a calculation 
with m max = 6 is shown in Fig. 1 for < L/N < 2. 







FIG. 3. Lines of constant density for L/N = 0.1,0.6,0.8, 
and 1.0 in a plane perpendicular to the z axis. On the darker 
curves |^| 2 = 0.3no, and on the lighter curves \^\ 2 = O.lno, 
where no = ff(z) 2 / 7ra osc- The unit of length is the oscillator 
length a osc . These pictures show how a vortex enters the cloud 
of bosons as the angular momentum per particle increases. 



We show the results with TO max = 6, since the occu- 
pancy of states with higher to is very low, and there- 
fore including such states would not alter the results on 
this scale. Figure 2 shows the corresponding interaction 
energy. Also Fig. 3 shows the lines of constant density, 
|*| 2 = constant for L/N = 0.1,0.6,0.8, and 1.0. Figure 3 
shows the gradual transition from mostly quadrupole and 
to a less extent octupole excitations, which are present 
at low angular momentum, to vortex-like structures as L 
approaches N. We should also mention that the struc- 
tures in Fig. 3, as well as those in Figs. 4, 5 and 6, rotate 
with an angular frequency fl given by 



_ 1 dE tot _ 1 d(V) 
h dL W Nh dl 



(14) 



which is lower than the trap frequency to. Here E tot is 
the total energy of the system. 

When L increases beyond N, the rotational invariance 
for L/N = 1 is lost. Density contours for various val- 
ues of L between N and 2N are shown in Fig. 4. These 
were calculated including the states up to to = 6. For 
L > 1.75JV the optimal wavefunction has a two- fold axis 
of symmetry, and the odd-TO coefficients in the wavefunc- 
tion vanish smoothly as the transition is approached, as 
shown in Fig. 1. In Figs. 5 and 6 we show contours for 
L = 2N and L — 2.17V. There is a first-order transi- 
tion from a state with two-fold symmetry to one with 
three- fold symmetry for L w 2.03/V. The solution for 
L/N — 2.0 with the states m = 0, 2, 4, 6, and 8 con- 
sidered has an energy of w 0.17577V 2 wo, whereas the one 
with the three-fold symmetry, with to = 0, 3, 6, and 9, 
has 0.1761/V 2 i;o. In contrast for L/N = 2.1, the state 
with the three-fold symmetry has an energy 0.1691iV 2 i;o, 
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FIG. 4. Lines of constant density for L/N = 1.25, 1.5, 1.75 
2.0, and 2.1 in a plane perpendicular to the z axis. On 
the darker curves |$| 2 = 0.3no, and on the lighter curves 
|*| 2 = O.lno. 



which is lower than the solution with the two-fold sym- 
metry, with an energy 0.170(W 2 vo- 

More generally, we have found that as L/N increases, 
the lowest-energy states are the ones where a vortex array 
is formed, in agreement with the results of Ref. Hj. 



where vq — Uo / |<I>o| 4 dr. One can easily see from 
Eq. ([!(]) that at this level of approximation the excita- 
tions with the highest gain in interaction energy per unit 
of angular momentum are the ones with A = 2 or A = 3, 
i.e., quadrupole or octupole excitations. 

We now calculate the interaction energy for low val- 
ues of The calculation of the energies of elementary 
excitations indicates that one would expect quadrupolar 
and octupolar modes to be the most important ones for 
small I. To determine the most energetically favorable 
way of giving the system angular momentum, one has to 
identify the behavior of | C2 1 2 and | C3 1 2 as I = L/N — ► 0, 
and then it is possible to build up a whole power-series 
expansion. Motivated by the fact that the A = 2 and 3 
excitations are degenerate, and are the ones which give 
the highest gain in energy per unit of angular momen- 
tum, we assume that both ci and c 3 are of order I 1 / 2 . 
As we show below, it is the mode-mode interaction that 
lifts this degeneracy, making the A = 2 mode dominant 
for low values of angular momentum. 

It is instructive to give an explicit example, so let us 
assume that we wish to examine the interaction energy up 
to order I 2 . In order to minimize the interaction energy 
to this order, the states with m = 1, 4, 5, and 6 need 
to be considered, since the phases of off-diagonal terms, 
like for example | cq | j Cx j | C2 1 1 C3 1 , can be chosen to have a 
negative sign, and thus lower the energy as compared to 
the case where only cq, C2, and c 3 are non-zero. A useful 
formula for the matrix elements of the potential is 



Jfc-t-Z,m+n 



(r)* m (r)* B (r)dr 
[k + l)\ 



2( fc +0V kill mini 



$o(r)| 4 dr. (17) 



B. Analytical results for L/N — > 



As will become clear below, to calculate the energy up to 
order I 2 we must include the following terms: 



We now turn to an analytic approach to the problem. 
One can systematically develop a power-series expansion 
for the occupancies |c m | 2 of the states, as well as for the 
energy in certain limits. We start with the case of very 
low angular momentum, I — L/N — ► 0. Working in terms 
of quantum-mechanical states, one of us has calculated 
in Ref. |j| the difference in the energies of two states 
where in the one state all N particles have m = and in 
the other state a particle is promoted to the state with 
m = A, and N — 1 particles have m = 0. If we denote 
the states as 



|0 J ' 



No ^JVi 2 N2 



(15) 



where N rn is the number of particles with angular mo- 
mentum mh, the two states are 10^), and lO^ -1 ^ 1 ). 
The difference e\ in the energy between these two states 
corresponding to this 2 A -pole excitation is given by 



€\ = Xtioj — 1 



1 



)A-1 



Nvo + O(vo), 



(16) 



y = (^M 4 + i| C ol 2 M 2 + i|c | 2 M 2 



+ 16 |C2 ' 



32 



C3 



C2 C 3 



IcoHcil 2 



V3, 



CO Ci C 2 C 3 



^|co| 2 |c 4 | 2 -^|co||c 2 | 2 |c 4 | 



+ ^|co| 2 |c 5 | 2 -^| C o||c 2 ||c 3 ||c 5 | 



|c ||c 3 | 2 |c 6 | ) N 2 v + O(Nv ). (18) 



In the above expression we have chosen the phases <j> m of 
the variational coefficients c m in such a way as to min- 
imize (V), and in the specific example we can arrange 
them so that all the off-diagonal matrix elements are 
negative. One of the phases can have any value, and 



4 



we make the choice 4>(s = 0. The rest of them can be 
expressed in terms of, say, <fti . We have found that up to 
to = 6 [to ^ 0] the expression 



4>m = ™/>l + (TO + 1)7T 



(19) 



gives the lowest energy. 

It is convenient to introduce the variable X = | C2 1 2 
| C3 1 2 , which is linear in I to leading order, and make 1 
of the constraints given by Eqs. (03) and (fL3j) to get 



I co I 



1 -X 



C4 



C 6 



|c 2 | 2 - 3A - J + | Cl | 2 + 4|c 4 | 2 + 5|c 5 | 2 + 6 1 c 6 1 2 ; 
|c 3 | 2 = I - 2X - | Cl | 2 - 4|c 4 | 2 - 5|c 5 | 2 - 6|c 6 | 2 . 

Then Eq. (|lj) takes the form 

w \ 2 4 16 8 32 



ci e 2 c 3 



(20) 




N 2 v + O{Nv ). 



(21) 



The last four terms in the above equation can lower the 
energy to order I 2 . For c 4 , for example, the energy is 
minimized if [see the third line of Eq. @] 



d 1 



9|c 4 | 



C4 



C4 



9 V6, |2| , 

C2 C4 , 



9|c 4 | 8 
|c 2 | 2 oc Z 



V6, ., 



(22) 



(23) 



Due to the non zero value of C4 the energy is lowered by 
an amount 



AS = 



1 



C4 - 



V6 



|c 2 | 2 |c 4 | JV 2 uo 



-^-\c 2 \ 4 N 2 v <xl 2 N 2 v . 
lb 



(24) 



It is remarkable that the term A£ exactly cancels the 
term 3 1 C4 1 4 / 1 6 in the second line of Eq. (|l8|). In a similar 
way Ci, C5, and c 6 can be expressed in terms of c 2 and c 3 
(and thus X), and Eq. ( [l8| ) takes the form of the effective 
Hamiltonian 



(V) 



34 



C3 



i| C 0| 4 + i| C 0| 2 |c 2 | 2 + i| C0 | 2 |c 3 | 2 



c 2 | 2 M 2 iV 2 W0 + O(7V U0 ), 



(25) 



or equivalently 



27 
17 



X 2 w 



-O{Nv ). (26) 



Minimizing the above expression with respect to X we 
find that X — 1/2 and thus the angular momentum has 
to be carried by the m = 2 state alone, since |c 2 | 2 = 1/2 
and I C3 1 2 = up to terms linear in I. Also the quadratic 
correction to (V) vanishes. Therefore for L/N — * 0, the 
quadrupole (A = 2) excitations are dominant. This is one 
of the important conclusions of our study. We show in the 
Appendix that a diagrammatic perturbative expansion 
which assumes that only the states with to = 0, 2, and 3 
are occupied by a macroscopic number of particles, while 
all the other states are not, [but still contribute to the 
energy] gives the same result. If one goes to higher order 
in / the interaction energy has within the perturbative 
scheme a term of the form | C2 1 3 1 C3 1 2 , which includes all 
the processes that convert three A = 2 excitations to two 
A = 3 excitations. This term can combine with the term 
I C3 1 4 , which implies that it is possible that | C3 1 2 oc |c 2 | 3 oc 
I 3 / 2 , which actually turns out to be the case. Then c\, for 
example, is given according to the second line of Eq. (J2l" 
by 



|ci|=V3|c 2 || C3 |oc/ 5 / 4 . 
Using similar arguments we find that 

|c m | 2 oc Z m/2 forTO ^ 1, and | Cl | 2 oc l 5/2 . 
The leading terms in |c m | 2 are given by 



(27) 



(28) 



and 



I CO 



ci 



c 2 



|c 3 | 



c 4 



|c 5 | 



|c 6 | 



1 



1 



1 



2 3 
l 5/2 + 2l 3 , 

-I -I 3 ' 2 , 
2 

2^3/2 

3 



I 3 ' 2 , 



1173 



- I 2 - - I 5 ' 2 
8 2 816 

-Z 5 / 2 --/ 3 , 
15 15 ' 

144 ' 



and the corresponding interaction energy is 



(V) 



1 I 



0(l 4 ) 



N 2 v + O(Nv ). 



(29) 



(30) 



The above equation is another basic result of our study, 
namely that the interaction energy drops linearly with 
the angular momentum up to the order we have exam- 
ined, for L/N — > 0, in agreement with our numerical 
simulations and with those of Refs. 
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C. Analytical results for / » 1 

We now turn to the region L /N « 1 . When the angular 
momentum per particle is exactly equal to 1, the lowest- 
energy state is the one where |ci| 2 = 1, and corresponds 
to a vortex state. We consider the two cases I < 1 and 
I > 1 separately, starting with I < 1. 

a. Analytical results for I < 1 

The simplest way to create a state with I < 1 from that 
with I = 1 is to transfer particles from the m = 1 state 
to the ground state. However, the energy can be even 
lower if also the m = 2 state is populated, because of the 
off-diagonal term |co||ci| 2 |c2|. The interaction energy up 
to order I = 1 — L /N is found by minimizing the potential 
energy, retaining only the coefficients cq, ci, and c%. This 
is 

W=Q| C i| 4 + |co| 2 | Cl | 2 + J| Cl | 2 |c 2 | 2 

-^\cq\\ Ci \ 2 \c 2 \\n 2 v + O{Nv ), (31) 

where we have used the fact that for this case too the 
phases may be shown to be given by Eq. (jl9|). Equation 
( |l9| ) is valid for small negative L/N — 1 at least up to 
m = 4. Thus in this limit 



and the interaction energy to order P is, 



|co| 2 oc |c 2 | 2 oc I. 



(32) 



To obtain the coefficients of proportionality it is conve- 
nient to use the following parametrization: 



|c | 2 = (l+a)l, 
M 2 =1- (l+2a)Z, 
|c 2 | 2 = al, 



(33) 



where a is a variational parameter. Minimizing the inter- 
action energy in Eq. ( |3l| ) with respect to a we find that 
a = 1. More generally using similar arguments we find 
that to leading order 



IcJ 2 oc D m -^ 



(34) 



and the explicit expressions for the coefficients are 

M 2 = 2f-^P, 



and led 



N 2 = l-3[+f P, 
\c 2 \ 2 =T-\p, 

|C3| 2 = ^, 

p 



12' 



(35) 



(V) = 



N 2 v„ + O(Nv ). 



(36) 



The above equation implies that also in the region I < 1 
the interaction energy varies linearly with the angular 
momentum to the order we have examined, which is also 
in agreement with the numerical simulations. The coef- 
ficient of the linear term is the same as the one we found 
for small values of the angular momentum. 

Equations (^) and ( |36| ) as well as the numerical results 
[see Fig. 2] strongly suggest that the interaction energy 
(V) drops linearly as a function of L/N in the whole re- 
gion < L/N < 1. The same result was derived by 
Butts and Rokhsar Q] numerically within the mean-field 
approximation. In Ref. Bertsch and Papenbrock per- 
formed exact diagonalizations of degenerate states of a 
given L and found that up to machine accuracy the en- 
ergy of the lowest state for a given L varies linearly with 
L in the range 2 < L < N , in agreement with our analytic 
expansions. This result is specific to the contact form of 
the effective interaction, and is probably connected to a 
hidden symmetry in our Hamiltonian H, as discussed by 
Pitaevskii and Rosch where the same Hamiltonian 
was considered in the context of breathing modes. 

b. Analytical results for I > 1 

We turn now to the case I > 1. We calculate 
the difference in energy between the states ll^), and 
ll^ -1 , (A + l) 1 ), by a method similar to that which for 
small I led to Eqs. ( ff5| ) and (|l6|). The energy e\ of this 
2 A -pole excitation with L = N + A, with A <C N is given 
by 

e x = \Tw - ~ (l - Nv + O(v ). (37) 

This formula implies that for I > 1, the single-particle 
excitations with the lowest energy per unit of angular 
momentum are those with A = 4 or 5, which means that 
the actual angular momentum carried by the particles 
is m = 5 or 6. In contrast to the low-angular momen- 
tum case, here both | C5 1 2 and |ce | 2 vary linearly with I, 
I C5 1 2 oc I Ce 1 2 oc I, where I = L/N — 1. In addition, in this 
regime we find that the energy has corrections of higher 
order than linear. Using similar arguments to those given 
before for small I, we find, to order I 2 , 



l c o| = 

M 2 = 
M 2 = 
M 2 = 
M 2 = 
M 2 = 



0.1213 P, 

1 - 0.2241 1, 

0.1934 P, 

0.1205 1, 

0.1036 1, 

1.7 x 10~ 3 P, 
3 72 



\c w \ = 1-3 x 10" J r, 
|cn| 2 = 8.6 x 10~ 3 P, 



(38) 
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and the interaction energy is 



(V) = 



I _ ^ + 6.7 xlO- 3 P + 0(P) 
4 64 



N 2 v + O(Nvq). 

(39) 



If we compare the expressions ([36|) and ( p9| ) for the inter- 
action energy we see that there is a change in its slope, 
d{V)/dL, as L passes N, from -7Vu /4 for L < N to 
-5Nv /6A for L > N. 



nodes in the density reflects the presence of two displaced 
vortices, and thus we see that the lowest-energy state of 
the system has two separated vortices and not a doubly 
quantized vortex. This clearly demonstrates the instabil- 
ity of the double-quantized vortex state to formation of 
two vortices plus surface waves. 



IV. BEYOND THE MEAN-FIELD 
APPROXIMATION 



D. Results for higher values of L/N 

We mentioned earlier that as the angular momen- 
tum per particle increases even further, there are certain 
ranges of values of L/N over which the state with the 
lowest energy has a specific symmetry. The lowest value 
of L/N for which this occurs is rs 1.75 and the symme- 
try of the state is two-fold, i.e., only ci m ^ 0. We have 
examined analytically as an example the case L/N = 2. 
Keeping only the first three non-zero coefficients, which 
are the dominant ones, we find to order I = L/N — 2 that 



and 



|c 2 | 

|C4| 



\co\ 2 =(3(l)-l/4, 



1 - 2/3(0, 
/3(f) + f/4, 



(40) 



where 



/3(f) 



3092-48V6 17(64^+109): 



12695 
0.2343 + 0.4449^, 



10156 



(41) 



or 



I co | 
leal 

\C4\ 



0.2343 + 0.1949Z, 
0.5314 - 0.8897Z, 
0.2343 + 0.6949f, 



and the corresponding interaction energy is: 

(V) = [A- Bl + CP + O (P)] N 2 v + O(Nv ), 
where 



(42) 
(43) 



A = 
B 



_3_ 

16 " 
1 

512 



/3o 
32 



(7 - 4V6) 



4 



^(64V6-109), 



(4 + 85/3 ), 



(44) 



and (3q = /3(0). The actual numbers which appear in 
Eq. ph are 



(V) « [0.1773 - 0.0467/ + 0.0170Z 2 + O (F)] N 2 v 

+O(Nv ). (45) 

There is no change in the slope of the interaction energy 
as L passes 2N. Figure 5 shows lines of constant den- 
sity, |^/| 2 = constant, for L/N = 2. The occurence of two 



Another way of approaching the problem of rotation, 
is to diagonalize the Hamiltonian within the space of de- 
generate states. This approach goes beyond the mean- 
field approximation, since in mean field theory the many- 
body wavefunction is the product of the single-particle 
states, whereas the diagonalization allows for the many- 
body state to have all kinds of correlations between the 
particles. This technique can be used by taking into ac- 
count the whole set of states JjJ , but it is convenient and 
pedagogical to work in a restricted space, appropriately 
chosen. 

As an example we consider small negative L/N — 1. 
From the analysis of Sec. Ill C we know that in this limit 
the states with m — 1,0, and 2 are dominant. Therefore 
the eigenstates 



|/i,f) = l^-'- 2 /, 2") 



(46) 



with N particles and L = N — I units of angular momen- 
tum are expected to provide a good basis for describing 
the low-lying states for / <§; N. We restrict ourselves to 
this limit and demonstrate how one can derive an effec- 
tive Hamiltonian which can be diagonalized exactly. In 
the limit we consider, I is <C N, and thus /i ~ I <C N. 
The diagonal matrix elements in the Hamiltonian are, up 
to terms of order N, 



±N(N-l) + ±lN+^N ] ,•„. 



(47) 



and the off-diagonal matrix elements are 



72- 



( M + l\V\n) « ^Nv V (/i + I + 1)(M + !)■ (48) 

Ignoring for the moment the (diagonal) first term of 
Eq. p?]), which corresponds to the interaction energy of 
the state 11^), we see from Eqs. ( f47| ) and fl4S| ) that we 
have to diagonalize the Hamiltonian 



H = 



1 t 1 t 

2 a a + ^ a 2 a 2 



V2, 



(<4aJ + a 2 a ) 



Nv , (49) 



which can be done exactly by use of a Bogoliubov trans- 
formation. Here a m is an annihilation operator that de- 
stroys a particle with angular momentum mh. Introduc- 
ing the operators c and d given by 
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c — a\ + V2a 2 ; d — V2a + a|, 

we may write the Hamiltonian as 

N + 
H = —(dU-l)v . 



(50) 



(51) 



Acting on states of the type (^) the operator Sd — c^c 
is diagonal, and has an eigenvalue N — L. Therefore <fd 
can be eliminated and from Eq. (EW we obtain 



H=*-(N-L-l)v + ^ 



N , 
1 



c'c. 



(52) 



The total interaction energy is thus the eigenenergy of H 
plus the diagonal part N(N — 1)uq/4, or 

(V) = j(N- l)v + j(N -L-1+ (c+c))«o + O(v ) 
= N ( 2N - L - 2 K + \{Jc)Nv + O(v ). (53) 

In the ground state (c'c) — 0, and the energy given 
by Eq. (p33j) is the same as that derived numerically in 



Ref. [7j. The presence of the term (c^c) in Eq. fl53| ) im- 
plies that the excited states in this limit of small nega- 
tive L /N — 1 are separated from the ground state by an 
amount 



N 

AE = —va + O(vo). 



(54) 



We have also performed numerical diagonalization, and 
we have confirmed the above result ([53]), as well as 
Eq. (|54|). The average occupancy of the states of the 
non-interacting problem is, for the lowest-energy state 
and L = N: 



1 



(55) 

and thus in the limit N — > oo there is agreement between 
the mean-field approximation and the present one. 



V. SUMMARY AND CONCLUSIONS 

To summarize, we have studied the lowest-energy 
states of a system of rotating, weakly interacting harmon- 
ically trapped bosons. Within the mean-field approxima- 
tion, for L/N — > we find that the angular momentum 
is carried mainly by quadrupole (|m| = 2) excitations. 
We have demonstrated that diagrammatic perturbation 
theory also leads to the same results as the method we 
have used here. 

For L/N = 1 the angular momentum is carried by 
particles in the m = 1 state, while for small negative 
L /N — 1 the to — and to = 2 states are also populated. 
In the limits L/N — ► and L/N — ► 1 the energy is a 



linear function of the angular momentum up to the order 
we have explored, while numerically this linearity persists 
in the whole region < L/N < 1. This result is specific 
to the contact form of the effective interaction, and does 
not hold for more general interactions. 

For small positive L/N — 1, the states which carry the 
additional angular momentum are those with to = 5 and 
to = 6. In addition, as L passes N the derivative of 
the interaction energy with respect to the angular mo- 
mentum changes abruptly. We have also found that for 
L/N w 1.75 there is a second order phase transition 
and for 1.75 < L/N < 2.03 the lowest energy state has 
two-fold symmetry. At L/N 2.03, there is a first or- 
der phase transition to a state with three-fold symmetry. 
More generally, for higher values of L/N a vortex array 
develops. 

The Gross-Pitaevskii wavefunction is a power series in 
z = x + iy. Thus if one truncates the series at to = m ma x, 
the wavefunction will have TO max nodes. In the vicinity 
of a node at z = zq the wavefunction varies as z — zq and 
therefore each node corresponds to a singly quantized 
vortex having the same sense as the total angular mo- 
mentum. It is instructive to study how the vortex lines 
move as the angular momentum is increased. For low an- 
gular momentum, the condensate wavefunction has only 
to = and m = 2 components, and it is therefore pro- 
portional to [l-(«V2/2)(p/a osc )2 e 2#] e W-p 2 /2a^ c ) for 
the choice of 4>2 = 7r, according to Eq. ( |l9| ) with <j)\ = 0. 
This has vortices on the x axis at x — ±(2/Z 1 / 2 ) 1 / 2 a OS c- 
With increasing angular momentum, components of the 
wavefunction with odd to grow, and the two-fold sym- 
metry of the cloud is broken, as may be seen in Fig. 1 
for L/N = 0.1, one of the vortices moving to larger dis- 
tances, and the other to smaller ones. The C3 term leads 
to a third vortex at large distances from the origin. For 
L/N = 1 there is only one vortex, which is at the ori- 
gin. With further increase in L/N, the velocity field is 
at first still dominated by a vortex close to the origin, 
but subsequently a second vortex moves into the cloud 
until at L/N m 1.75 the two-fold symmetry is restored. 
As L/N increases towards the value 2.03, at which the 
first-order transition to the state with three-fold symme- 
try mentioned above occurs, the separation of the two 
vortices changes little. 

In this paper we have also investigated effects not in- 
cluded in mean-field theory by diagonalizing a model 
Hamiltonian for L close to, but less from, N. We find 
that for L = N the occupancy of the to = 1 state is 1, 
with corrections of order 1/N. We have calculated the 
energy up to terms of order N. Finally we also found 
that the low-lying excited states are separated from the 
lowest state by energies of order A^^o of the same angular 
momentum. 

In this study we have examined the limit of weak in- 
teractions. When the interaction energy per particle nUo 
becomes comparable to or greater than hu>, components 
of the wavefunction that are not members of the low- 
est multiplet in the absence of interaction must be in- 
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eluded. Calculations for this regime based on the Gross- 
Pitaevskii equation have been carried out by Isoshima 
and Machida (l| . Comparison of our results with theirs 
is difficult because these authors calculated the lowest 
energy state in a rotating frame, rather than the lowest 
energy state for a given angular momentum. 

One question of importance both conceptually and be- 
cause of its relevance to experiment is whether or not 
the states are stable to small perturbations, and if they 
are not, what is the lifetime of the state. The answer 
to these questions depends on the nature of the pertur- 
bation, whether it is due to a deformation of the trap, 
or to interactions with particles outside the condensate, 
and we shall discuss it elsewhere. 

In our calculations above we have shown for a particu- 
lar example that the Gross-Pitaevskii approach gives cor- 
rectly the contribution to the energy of order iV 2 . This 
result, which is alluded to in Ref. [Q, is more general, 
and in a future publication |l4| it will be shown how the 
Gross-Pitaevskii approach is recovered as the first term 
in an expansion in powers of 1/N. The method may be 
extended to calculate contributions to the energy of order 
TV which are in excellent agreement with results obtained 
by numerical diagonalization of the Hamiltonian. 



APPENDIX: PERTURBATION THEORY 
APPROACH 



which is precisely the correction A£ given by Eq. ( |24| ) 
(plus terms of order Nvq). Similarly the other dia- 
grams shown below give — 5|c3| 4 /544, — | C2 1 2 1 C3 1 3 / 8 , and 
— 3 1 C2 1 2 1 C3 1 2 /4 in units of N 2 vo, respectively, and are 
identical to the corrections given by the terms in the last 
four lines of Eq. (|l^ ) . 



m=2 m=2 m=3 m=3 m— 2 m— 3 m=2 



m=3 







m=2 m=2 m=3 m=3 m=2 m— 3 



m=3 



m=2 



FIG. 5. The four diagrams contributing to the interaction 
energy to order I 2 . 
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We show in this appendix that one can use pertur- 
bation theory to derive an effective Hamiltonian in the 
region L/N — > 0, corresponding to Eq. (|25|). We assume 
that only the states with m = 0, 2, and 3 are macroscop- 
ically occupied. However, other states (the m = 1,4,5 
and 6 in this case) give corrections to the energy that 
can be treated perturbatively. Let us demonstrate how 
this works by considering the interaction energy up to 
I 2 . As long as both C2 and C3 vary as Z 1 / 2 , the only pro- 
cesses that contribute to the interaction energy up to I 2 
are shown in Fig. 7. 

Let us consider the first process on the left as an exam- 
ple. The matrix element M corresponding to the vertex, 
where two particles with m = 2 scatter to states with 
m = and m — 4, is equal to 



16 



(Al) 



From Eq. (|l6| ) we find that the difference in the energy 
between the intermediate state and the initial state is 



Se 



Nv 



+ O(vo), 



(A2) 



and thus perturbation theory implies that the correction 
to the energy is 

= -^N 2 v = -^| C2 | 4 A^o, (A3) 
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